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ABSTRACT 

This is the first in a series of papers in which we study the application of spectroastrometry in the context of gas kinematical studies 
aimed at measuring the mass of supermassive black holes. The spectroastrometrical method consists in measuring the photocenter of 
light emission in different wavelength or velocity channels. In particular we explore the potential of spectroastrometry of gas emission 
lines in galaxy nuclei to constrain the kinematics of rotating gas disks and to measure the mass of putative supermassive black holes. 
By means of detailed simulations and test cases, we show that the fundamental advantage of spectroastrometry is that it can provide 
information on the gravitational potential of a galaxy on scales significantly smaller (~ 1/10) than the limit imposed by the spatial 
resolution of the observations. We then describe a simple method to infer detailed kinematical informations from spectroastrometry 
in longslit spectra and to measure the mass of nuclear mass concentrations. Such method can be applied straightforwardly to integral 
field spectra, which do not have the complexities due to a partial spatial covering of the source in the case of longslit spectra. 

Key words. Line: profiles - Techniques: high angular resolution - Techniques: spectroscopic - Galaxies: active - Galaxies: kinemat- 
ics and dynamics - Galaxies:nuclei 



1. Introduction 

One of the fundamental open questions of modern astrophysics 
is understanding the physical processes that transformed the 
nearly homogeneous primordial medium into the present-day 
universe, characterized by a wealth of complex structures such 
as galaxies and clusters of galaxies. Understanding how galax- 
ies formed and how they become the complex systems we ob- 
serve today is therefore a major theoretical and observational 
effort. In the last few years, strong evidence has emerged for 
the existence of tight links between supermassive black holes 
(BH), nuclear activity and galaxy evolution. These links reveal 
what is now called as co-evolution of black holes and their 
host galaxies. Strong evidence is provided by the discovery of 
'relic' BHs in the center of most nearby galaxies, and that BH 
masses (Mbh - 10 6 - 10 10 M ) are tightly proportional to struc- 
tural parameters of the host spheroid like mass, luminosity and 
stellar velocity dispersion (e.g. Kormendy & Richstone 1995, 
Gebhardt et al. 2000, Ferrarese & Merritt 2000, Marconi & Hunt 
2003, Ferrarese & Ford 2005, Graham & Driver 2008 and refer- 
ences therein). Moreover, while it is widely accepted that Active 
Galactic Nuclei (AGN) are powered by accretion of matter on a 
supermassive BH, it has been possible to show that BH growth is 
mostly due to accretion of matter during AGN activity, and there- 
fore that most galaxies went through a phase of strong nuclear 
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activity (Soltan 1982, Yu & Tremaine, Marconi et al. 2004). It 
is believed that the physical mechanism responsible for this co- 
evolution of BHs an their galaxies is probably the feedback by 
the AGN, i.e. the accreting BH, on the host galaxy (Silk & Rees 
1998, Fabian 1999, Granato et al. 2004, Di Matteo et al. 2005, 
Menci et al. 2006, Bower et al. 2006). 

In order to proceed further it is important to secure the most 
evident sign of co-evolution, the correlations between BH mass 
and galaxy properties which can be achieved by increasing the 
number, accuracy and mass range of existing measurements. 
Supermassive BHs are detected and their masses measured by 
studying the kinematics of gas or stars in galaxy nuclei and, cur- 
rently, about ~ 50 BH mass measurements most of which in 
the 10 7 - 10 9 M Q range and only very few measurements below 
and above those limits (see, e.g., the most recent compilation by 
Graham 2008). 

This paper, the first in a series, deals with BH mass measure- 
ments from gas kinematics and presents a new method, based 
on spectroastrometry, which can provide a simple but accurate 
way to estimate BH mass and which partly overcomes the limi- 
tations due to spatial resolution which plague the 'classical' gas 
(or stellar) kinematical methods. 

In Sect.|2]we introduce the standard method for gas kinemat- 
ical studies, that is based on the gas rotation curves, and briefly 
discuss its characteristics and limitations. In Sect. [3] we intro- 
duce the new gas kinematical method based on spectroastrom- 
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etry. We explain the basis of this approach that consist on mea- 
suring "spectroastrometric curves" (Sect. |3.1| l and show simu- 
lations based on a model of the gas dynamics from which we 
want to learn how the spectroastrometric curve changes in func- 
tion of some parameters of the model (Sect. 3.2 1. In Sect. |4]we 
explain the practical application of the method. We start by pre- 
senting a method for using simultaneously several spectroastro- 
metric curves of the same source (Sect. |4.1| i and we consider 
its application on noisy data (Sect. |4.2| i. In Sect. [5] we present 
a trivial fitting method for recovering the values of the various 
model parameters using long slit spectroscopy with noisy data. 
Finally, in Sect. [6] we describe the practical application of the 
method with Integral Field Units (IFU's) and in Sect. [H] we draw 
our conclusions. In appendix[A]we discuss in details of the spec- 
troastrometric measurements with special regard to the determi- 
nation of the light centroids and in appendix [B] we describe in 
detail the method to recover 2D spectroastrometric maps from 
multiple slit spectra. 

2. The standard gas kinematical method 

The standard method of gas kinematics is based on spatially re- 
solved spectra of the nuclear region of a galaxy. This method 
consists in recovering the rotation curve of a given gas emission 
line from a longslit spectrum. The modeling assumes that the gas 
is rotating in a thin disc configuration (neglecting hydrodynam- 
ical effects) under the effect of the gravitational potential of the 
stellar mass and of a pointlike dark mass Mbh that is the BH. 
The value of Mbh and other unknown parameters of the model 
are obtained by a fitting the rotation curves (see |Marconi et al. 
(2006) and references therein for details). Because of the many 
unknown parameters of the model, to better constrain the fit, usu- 
ally many spectra of the galactic nucleus are obtained each with 
a different orientation of the spectroscopic slit and a simultane- 
ous fit of all of them is performed. 

Clearly, the ability to detect the presence of a BH and mea- 
suring its mass strongly depends on the signal-to-noise ratio of 
the data and only poor constraints on Mbh can be obtained from 
low S IN data. However, even with high S /N data, the fundamen- 
tal limit of the standard gas kinematical method ("rotation curves 
method" hereafter) resides in the ability to spatially resolve the 
region where the gravitational potential of the BH dominates 
with respect to the contribution of the stars. Similar consider- 
ations apply equally to stellar dynamical measurements. 

As a rule of thumb, this region corresponds to the so-called 
sphere of influence of the BH. The radius of the sphere of influ- 
ence (rsH) can be estimated as (Binney J.|1987|): 
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where G is the gravitational constant, Mbh the BH mass and 
cr* is the velocity dispersion of the stars in the galaxy. For a 
galaxy with BH mass Mbh ~ 1O 8 M0 and stellar velocity dis- 
persion (T star ~ 200km /s we obtain rBH ~ H-2pc which, for a 
very nearby galaxy at distance D ~ 3 Mpc, provides an apparent 
size of ~ 0.7" . If a typical very good spatial resolution available 
from ground based observations is of the order of ~ 0.5" (this 
Full Width Half Maximum (hereafter FWHM) value should then 
be compared with twice rBH), we can notice that the sphere of in- 
fluence is only marginally resolved even for very nearby galaxies 
with moderately large BHs. For a galaxy distance of 30 Mpc the 
sphere of influence becomes marginally resolved even with the 
Hubble Space Telescope (HST) which provides the best spatial 



resolution currently available from space. Therefore, the rota- 
tion curves method can detect only BHs with moderately high 
masses and located in nearby galaxies; this is a strong limita- 
tion from the point of view of a "demographic" study of BHs in 
galactic nuclei, because the adopted "investigation tool" cannot 
reveal the entire BH population. 



3. The spectroastrometric gas kinematical method 

The use of sp e ctroastrometry was origi nally introduced by 
Beckers| ([1982]), |Christy et al!] ( fT983] ) and |Aime et al.| ( fT988> 
to detect unresolved binaries. However these earlier studies re- 



quired specialist instrumentation, and it was not until the work 
of Bailey (1998) that this method was exploited using standard 
common user instrumentation: a longslit CCD spectrograph. 

Subsequently spectroastrometry has been used by several au- 
thors to study pre main sequence binaries (Bai nes et al.| [2004 



Porter e t al. 2004] [2005] ) and the presence of inflow or outflow 



or the disk structure on the gas surrounding pre main sequence 
stars (|Takami et al.||2003| |Whelan et aL]|2005| ). More recently, 
Brannigan et al. (2006 ) discussed the presence and the detection 
of artifact in the output of the method and Pontoppidan et al. 



(2008) used spectroastrometry on 4.7 pmn CO lines to study the 
kinematical properties of proto-planetary disks. 

The fundamental advantage of the spectroastrometric 
method is that, in principle, it can provide position measure- 
ments on scales smaller than the spatial resolution of the ob- 
servations. 

We will now explain the general principle of the spectroas- 
trometric method with a simple example: consider two point-like 
sources located at a distance smaller than the spatial resolution 
of the telescope; these sources will be seen as spatially unre- 
solved with their relative distance not measurable from a conven- 
tional image. However, if in the two sources are present spectral 
features, such as absorption or emission lines at different wave- 
lengths, the light profiles extracted from a longslit spectrum at 
these wavelengths will show the two sources separately. From 
the difference in the centroid of the light profiles at these two 
wavelengths one can estimate the separation between the two 
sources even if this is much smaller that the spatial resolution. 
This "overcoming" of the spatial resolution limit is made possi- 
ble by the "spectral" separation of the two sources. 

While this method has been applied to unresolved pointlike 
sources as binary stars and protostellar systems, it has never been 
applied to the problem of measuring BH masses from the nu- 
clear gas emission in galaxies. Motivated by the need to over- 
come the shortcomings of the standard gas kinematic method 
described above in Sect. [2] in this paper we study the applica- 
tion of the spectroastrometric method to gas kinematical studies 
of the mass of BH's in galactic nuclei. We will start presenting 
the spectroastrometric method based on the same longslit spec- 
tra used for the rotation curves method but instead of the rotation 
curves we recover the "spectroastrometric curve". However, as 
shown in Sect. [6] the application of spectroastrometry to inte- 
gral field data is even more simple and straightforward than to 
longslit spectra. 

Here we focus on the theory of the method and the devel- 
opment of a practical framework for its application, exploring 
its capabilities and limitations using simulated data in order to 
understand how the spectroastrometric curve is affected by the 
object itself or by the instrumental setup. In subsequent papers 
we will apply the method to the kinematical data of real galaxies 
to yield new improved BH mass measurements. 
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Fig. 1. Example of the measurements on which the construction of rotation and spectroastrometric curves is based. Upper left panel: 
the longslit spectrum of a given continuum-subtracted emission line (isophotes) with superimposed two particular sections at fixed 
slit position and velocity (dotted lines). Bottom left panel: the profile relative to the Sect, at fixed slit position with superimposed 
the position of its centroid (dashed line). Upper right panel: the profile relative to the Sect, at fixed velocity with superimposed 
the position of its centroid (dashed line). Bottom right panel: the position-velocity diagram (PVD) of the upper left panel with 
superimposed rotation (blue dots with error bars) and spectroastrometric curves (red points with error bars). In all panels velocities 
are measured relative to the systemic velocity of the galaxy Vsys which is set to the zero point of the velocity scale. 



In this paper we concentrate on the application of the spec- 
troastrometric method to continuum subtracted spectra in or- 
der to focus exclusively on the gas kinematics. Indeed, among 
other things, the underlying continuum only dilutes or modifies 
the spectroastrometric signal expected for the spatial distribution 
and kinematics of emission line gas. 

3. 1 . The spectroastrometric curve 

Here we explain how the spectroastrometric curve is obtained 
and its main differences with "classical" rotation curves. 

A longslit spectrum of a continuum subtracted emission 
line provides a pixel array whose axes map the dispersion and 
slit directions, the so-called position-velocity diagram (here- 
after PVD). The slit axis maps the observed position along the 
slit. The dispersion axis maps the wavelength of emission from 
which the line of sight velocity can be derived. The upper left 
panel in Fig. [T] shows the isophotes of an emission line from 
a PVD. Ideally, for infinite S/N and perfectly circularly rotat- 
ing gas, the rotation curve denotes the mean gas velocity as a 



function of the position along the slit. In practice, that curve is 
obtained by fitting the observed line profiles along the slit with 
gaussian functions which provide an estimate of the average ve- 
locity after discarding components which are clearly not circu- 
larly rotating. The lower left panel of Fig. [T] displays the line 
profile extracted at the slit position marked by the horizontal dot- 
ted line in the PVD. The spectroastrometric curve provides the 
mean position of the emitting gas as a function of velocity. In 
ideal data that curve is derived by taking the light profile of the 
line at given velocities along the dispersion direction and mea- 
suring the corresponding mean emission centroids. In practice, 
the finite S/N of real data requires more complex measurements, 
which are described below. The upper right panel of Fig. |T|dis- 
plays the light profile extracted at the velocity marked by the 
vertical dotted line in the PVD. The two curves can be compared 
in lower right panel superimposed on the isophotes in the PVD 
diagram. Clearly, the two methods analyze the same spectrum 
from complementary points of view. 



4 



Gnerucci et aL: Spectroastrometry of rotating gas disks and supermassive black holes. 




-0.1 



o 
On 



Velocity (km/s) 




M 

O 



O 

0, 




0.3 



0.2 



r o.i 



r„ -o.o 



-o.i 



-0.2 



-0.3 



Velocity (km/s) 



-500 





500 






Velocity (km/s) 






■ ■ i ■ 


/ 1 *\ 1 \ 1 1 


■ i 

H 




i 


. \ i . . J 


1 pixel : 

i 


-500 





500 






Velocity (km/s) 







Fig. 2. Spectroastrometric curve (black solid line) obtained from a simulated spectrum for a rotating gas disk with various values of 
the central BH and stellar mass. In gray we plot the isophotes of the simulated line spectrum. Upper left panel: case of Mbh = 1O 8 M0 
BH without mass contribution from stars. Upper right panel: stars only, with Mbh = 0. Lower panels: Mbh - 1O 8 M0 (left) and 
Mbh = 10 6 ' 5 Mq (right) with included the same stellar mass distribution as in model in the upper right panel. Note that the differences 
in the kinematics observed in the left panels (both BH dominated, but with and without stellar contributions) are marginal, and only 
visible i.e. beyond 072 from the nucleus. 



3.2. Simulations 

In the following we show the results of the tests based on sim- 
ulations of rotating gas disks, outlining the effect of varying the 
free parameters of the simulation and of the spectroscopic obser- 
vations. 

The model we use in our simulations is based on a thin 
gas disk rotating in a plane under the effect of the gravitational 
potential produced by a mass distribution of stars in a galaxy 
nucleus and by a pointlike dark mass, the central BH Mbh- 
Therefore, we assume that the gas is circularly rotating in the 
disk plane with the rotational velocity uniquely determined by 
the combination of Mbh and the stellar mass distribution. This 
model depends on several parameters: the dynamical parame- 
ters (the BH mass, the shape of the mass density function of 
the stars and the systemic velocity of the galaxy) and the geo- 
metrical parameters that establish the position and orientation of 
the gas disk (the distance of the galaxy or the angular distance 
scale, the inclination of the disk plane with respect to the line 
of sight and the orientation of the line of modes of the disk). 
The model takes into account the effect of the shape of the in- 
trinsic light distribution of the emission line on the sky plane, 



that is modeled analytically with a combination of gaussian or 
exponential functions. The model takes also into account the ef- 
fect of the instrumental Point Spread Function that is modeled 
with a gaussian function with a given FWHM, and the effect of 
the others instrumental setup parameters: the spectral resolution, 
the slit width and position angle, the detector's pixel size. With 
this model we can simulate a longslit spectrum of a gas emission 
line in a particular galaxy nucleus (see Mar coni et al. 12006) and 
references therein for a detailed description). 

For simplicity, in order to illustrate the technique, we make 
use of a basic reference model which is meant to closely match 
the physical parameters of a BH of mass 10 8 Mq in an ellipti- 
cal galaxy in the local universe observed with a typical longslit 
spectrograph, like ISAAC at the VLT ( |Moorwood et al.p 999). 
The key model parameters are: 

- The distance of the galaxy is set to 3.5 Mpc that corresponds 
to an angular distance scale of 17 pc/ "(e.g. the same of the 
galaxy Centaurus-A). 

- The disk inclination is set to 35°. 

- The disk line of nodes position angle (with respect to the 
North direction) is set to 0° 
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Fig. 3. Left panel: spectroastrometric curves from spectra differing only in spatial resolution. Solid line: spatial resolution of 0.1". 
Dotted line: spatial resolution of 0.5". Dashed line: spatial resolution of 1.0". The vertical long-dashed lines denote the "high 
velocities". Right panel: corresponding rotation curves. 
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Fig. 4. Left panel: spectroastrometric curves for spectra differing in the value BH mass. Solid line: Mbh = 0. Dotted line: Mbh - 
10 6 5 Mq. Right panel: rotation curves for the same models. 



- The amplitude of rotational velocity in the disk due to the 
stellar mass component at ±1" is set to + 200 km/s. 

- The FWHM of the spatial PSF is set to 0.5" as for typical 
high quality ground based observations. 

- The spectral resolution is set to 10 km/s. 

- The detector's angular pixel size is set to 0.15" x 0.15", re- 
sulting in a spatial oversampling of ~ 3.5. 

- The spectrograph slit is usually set parallel to the disk line of 
nodes unless otherwise specified. 



We have chosen a basic set of parameters which results in a 
resolved BH sphere of influence with the aim of clearly showing 
the spectroastrometric features of disk rotation. In the rest of the 
paper we will of course consider more extreme sets of parameter 
values which will result in non-resolved BH spheres of influence 
in order to show the full power of spectroastrometry. 



3.2.1. BH mass 

We first consider the effect of changing the value of the central 
BH mass on the spectroastrometric curve. 

We simulated the spectrum of a particular rotating gas disk 
model with different values of the mass of the central BH, and 
then we derived the spectroastrometric curves from these simu- 
lated spectra (see Fig. [2]). 

The spectroastrometric curve in Fig.|2]is centered at the sys- 
temic velocity of the galaxy. In these noiseless simulations, used 
for illustration purposes, the spectroastrometric curves are drawn 
only when the mean flux of the light profile along the slit is larger 
than 10~ 3 of the maximum flux of the spectrum. This is meant 
to give a representation of the low signal-to-noise regions that 
should be used in real data for the measurement of the centroid 
of the light profile, without considering unrealistically low flux 
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Fig. 5. Left panel: simulated spectroastrometric curves from spectra differing only in slit position angle (referred to the disk line 
of nodes 6 s ut - Own)- Right panel: measured amplitude of the spectroastrometric curve (distance between the maximum and the 
minimum) as a function of the slit position angle referred to the disk line of nodes; the dashed line denotes a cosine function 
opportunely rescaled to match the value of the amplitude of 9 s n, - Olon = 0. 



levels. We defer a discussion of the practical cut off fluxes im- 



posed by the signal-to-noise in real data to Sect. 4.2 



In the Mbh - case (Fig. [2] upper right panel) the gas kine- 
matics is due only to the gravitational potential of the stellar 
mass distribution and the spectroastrometric curve is monoton- 
ically decreasing at increasing velocity. In all other cases (Fig. 
|2]other panels), the points at high velocities (v < -300km/s and 
v >300km/s) tend to approach the 0" value for increasing veloci- 
ties (we consider "increasing" referring to |v— V 5 y. s |; V sys = for 
this simulation), giving the curve the characteristic "S"-shape. 
The reason for this behavior is that when the gas kinematics is 
dominated by the BH's gravitational potential, the gas at high 
velocity is located close to the BH and the light centroid of the 
high velocity gas has to be near to the BH position. 

Operatively the "high velocity" points are those where the 
line emission along the slit is spatially unresolved. In Appendix 
[A] we describe how to find these points as a by product of cen- 
troid determination. We can also recall that the presence of a 
turnover in the spectroastrometric curve is the signature of the 
presence of a BH. Therefore, the "high velocity" range is also 
characterized by a ~ 75% drop of the spectroastrometric signal 
with respect to the maximum shift reached at the turnover point. 

The velocity in a circular orbit of radius r around a point- 
like mass (the BH) is well approximated by a keplerian law with 
v = (GMsH/r) 1 ^ 2 in the innermost region of the gas disk where 
we can neglect the contribution of the gravitational potential 
of the stellar mass distribution. Combined with the above con- 
siderations, this also suggests that the spectroastrometric curve 
asymptotic value at high velocities provides an estimate of the 
BH position along the slit. 

For increasing Mbh values the spectroastrometric curve ex- 
tends to higher velocities (higher values of |v - V sv .sl) since this 
has the effect of increasing the amount of emission in the high 
velocity bins. However the limited extension in the velocity axis 
of the spectroastrometric curve with real data is due to the pres- 
ence of noise. At velocities where the flux of the line is too low 
with respect to noise we will not be able to calculate a reliable 
value of the centroid of the light profile along the slit. 



In conclusion the spectroastrometric curve reveals the pres- 
ence of a pointlike mass contribution to the gravitational poten- 
tial when it shows an "S-shaped" structure, with a turn-over of 
the high velocities components that get closer to the 0". It can 
then be concluded that the information about the BH resides pre- 
dominantly in the "high velocity" part of the curve. 

3.2.2. Spatial resolution 

Here we consider the effects of the spatial resolution on the 
spectroastrometric curve. The spatial resolution is the width 
(FWHM) of the point spread function which, in the model, is 
approximated by a gaussian function. 

The results of this test are shown in Fig. [3] where we com- 
pare spectroastrometric and normal rotation curves. In the left 
panel we present the spectroastrometric curves for models dif- 
fering only in spatial resolution while in the right panel we dis- 
play the corresponding rotation curves. The spectroastrometric 
curves differ in the "low velocity" range (-300 km/s < v < 300 
km/s), as they show a steeper gradient and a large amplitude. 

However, in the "high velocities" range differences are neg- 
ligible, ~ 0.1 pix at most. In section 3.2.1 we showed how the in- 



formation on the BH mass is encoded in the high velocity range 
and, in particular, the presence of the BH is revealed by the fact 
that the spectroastrometric curve approaches 0" at high veloci- 
ties. The spectroastrometric curve at the "high velocities" is al- 
most unchanged by worsening the spectral resolution, leaving 
the BH signature unaltered. 

For the "standard" rotation curves the information about the 
BH is also encoded in the presence of points at high velocities 
at low distances from the center. However, a lower spatial res- 
olution the BH signature is effectively canceled in the rotation 
curves when the sphere of influence of the BH is not resolved. In 
this simulation the apparent radius of the sphere of influence is 
~ 0.6", which is resolved in the case of 0.1" spatial resolution, 
partially resolved for 0.5" case, and unresolved for 1.0". 

Clearly, even in the case of the spectroastrometric curves a 
better spatial resolution is desirable. In fact, a poorer spatial res- 
olution results in a broadening of the light profile along the slit 
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Fig. 6. Left panel: iso-velocity contour map of the line of sight velocity field on the plane of sky for a rotating gas disk with the of 
line of nodes along the x axis, with superposed a 0.5" wide slit. Right panel: comparison of spectroastrometric curves for simulated 
spectra differing only in slit width (we remind that we assumed a spatial resolution of (X'5). Solid line: 0.2" wide slit. Dotted line: 
0.5" wide slit. Dashed line: 1.0" wide slit. 



and in a decreased accuracy with which the centroid position 
(and consequently the BH mass) can be measured. Note, how- 
ever, that the spatial resolution is not a free parameter, since it 
can be measured directly from the data. Its effects can be mod- 
eled and taken into proper account. 

We have just shown that the the rotation curve changes dras- 
tically when varying the spatial resolution while the spectroas- 
trometric curve does not. Now we conclude this section with a 
last set of simulations to showing in a qualitative but more ac- 
curate way how this method can really allow us to overcome the 
spatial resolution limit. 

The results of these simulations are shown in Fig. |4] In the 
left panel we display the spectroastrometric curves for two mod- 
els with the same spatial resolution (0.5") but different BH mass; 
in the right panel we display the corresponding rotation curves. 

By considering only the rotation curves, the case Mbh - 
is effectively indistinguishable from that with M BH = 10 6 5 Mq 
and indeed the sphere of influence is unresolved in this case (the 
apparent radius of the sphere of influence is r B H ~ 0.02"). 

Instead, for this particular set of simulations, a M BH = 
10 6 5 Mq mass appears to be still distinguishable from the Mbh = 
case; in the high velocity range the differences between the two 
models reach a value of ~ 072, i.e. ~ 1.3 pixels. Such difference 
can be measured in real data, assuming that we can achieve a rea- 
sonable accuracy of ~ 1 pixel in the measure of the photocenter 
at these velocities. In conclusion, according to these simulations 
we are able to detect a BH whose apparent size of the sphere of 
influence is less than ~ 1/10 of the spatial resolution. 

3.2.3. Slit position angle 

Here we consider how the spectroastrometric curve is affected by 
changing the position angle of the slit. We simulated the spec- 
trum of a particular rotating gas disk model with different values 
of the position angle of the slit and then we get the spectroastro- 
metric curves from this simulated spectra. 

In the left panel of Fig. [5] we can see the comparison of the 
spectroastrometric curves relative to different values of the slit 
position angle referred to the disk line of nodes position angle 



{0 s ut - Olon where 6 s u t and Olon are the position angles of the 
slit and disk line of nodes respectively). 

We can see that for 9 s u t — Olon = (slit aligned with the disk 
line of nodes) the spectroastrometric curve has the maximum 
amplitude. For increasing values of s u t — Olon the amplitude de- 
creases reaching a null value for s u, - Olon = 90°. For values 
larger than 90° the curve invert itself. This is clearly a geomet- 
rical projection effect just like the one affecting the amplitude 
of normal rotation curves as a function of s n, - Olon- In the 
right panel of Fig.[5]we show the curve amplitude as a function 
°f 9 slit ~ ^Lcwwhich is well approximated by a cosine function 
(dashed line). 

3.2.4. Slit width 

Another parameter which influences the spectroastrometric 
curve is the width of the slit with which the spectra are obtained. 
Each point of the spectroastrometric curve represents the cen- 
troid of the light profile along the slit at a given velocity but one 
can only sample the fraction of light emitted by the gas at that 
specific velocity that is intercepted by the slit. 

In Fig.[6j left panel, we show the map of line of sight velocity 
field of a rotating gas disk. The gas in a given velocity bin lies in 
the locus delimited by two subsequent isovelocity contours. The 
impact of superimposing the slit on the iso-velocity contour map 
is to artificially truncate the spatial regions contributing to line 
emission at a given velocity that are not confined within the slit 
extension. This results in distorted and displaced photo-centers 
in velocity space. 

In Fig.|6| right panel, we show the spectroastrometric curves 
obtained from models differing only in slit width. As long as the 
slit width is smaller or equal to the spatial resolution of the ob- 
servations (which is, in this case, 075 FWHM) then differences 
among spectroastrometric curves in the high velocity range are 
negligible (e.g., slits with 0.2", 0.5" width). On the contrary the 
spectroastrometric curve is significantly affected for larger slit 
widths, even in the high velocity range, because of the inclu- 
sion of more extended emission. This comparison indicates that 
one should select for the observations a slit width smaller or at 
most equal to the spatial resolution of the observations. As al- 
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ready noted discussing the effects of spatial resolution, the small 
residual differences in the spectroscopic curves corresponding to 
various slit widths can be effectively modeled out. 

3.2.5. Spectral resolution 



0.3 F"~~ r 
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however derived trading-off with the level of signal-to-noise nec- 
essary to build well defined spectroastreometric curves. 

Finally, we note that the value of cr can be estimated by mod- 
eling the classical rotation curves. It is then possible to establish, 
a posteriori, whether the BH measurement is affected by such an 
effect and eventually to validate its value. 
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Fig. 7. Comparison of spectroastrometric curves for spectra dif- 
fering only in spectral resolution. Solid line: spectral resolution 
of 10 km/s. Dotted line: spectral resolution of 50 km/s. Dashed 
line: spectral resolution of 100 km/s. Dot-dashed line: spectral 
resolution of 150 km/s. 

We now focus on the effects of the finite spectral resolution 
of the observations on the spectroastrometric curves. In Fig. [7J 
we show models differing only in spectral resolution. By de- 
creasing the spectral resolution, the amplitude of the spectroas- 
trometric curves is decreased and they show, at a given spatial 
offset, higher velocities. Effectively, due to the spectral convolu- 
tion, the curves are stretched along the velocity scale. To under- 
stand this effect we can approximate the true velocity profile and 
the instrumental line profile with gaussian functions. The result- 
ing line profile is the convolution of these two functions and is 
therefore a gaussian function with standard deviation 



-4 



cr 2 + crl 



(2) 



where cr is the resulting standard deviation, cr is the standard 
deviation of the velocity profile and cr is the standard deviation 
of the instrumental response function. 

We here consider as a source of the line broadening the un- 
resolved rotation of the gas that originates from the fact that one 
is observing with finite spatial resolution and is then not able to 
spatially resolve the high velocity regions close to the BH. In our 
simulations, the line width due to unresolved rotation reaches 
cr ~ 250 km/s at the galaxy's center. 

As expected, the spectroastrometric curves are almost un- 
changed at high spectral resolution (cr = 10 , 50 km/s) but they 
are significantly altered for the lower spectral resolutions consid- 
ered (cro = 100 , 150 km/s) when cr approaches the value of cr. 
As a consequence, the BH mass estimates derived from data of 
insufficient spectral resolution are systematically overestimated, 
since at a given position, one measures an artificially increased 
velocity. 

This result underscores the importance of using data of as 
high as possible spectral resolution. The optimal value must be 



Fig. 8. Comparison of spectroastrometric curves for spectra dif- 
fering only in intrinsic velocity dispersion. Solid line: intrinsic 
velocity dispersion of 10 km/s. Dotted line: intrinsic velocity 
dispersion of 50 km/s. Dashed line: s intrinsic velocity disper- 
sion of 100 km/s. Dot-dashed line: intrinsic velocity dispersion 
of 150 km/s. 



3.2.6. Intrinsic velocity dispersion 

Another source of line broadening might be due to non circular 
or chaotic motions in the gas which can be modeled by adding 
an "intrinsic" velocity dispersion to the gas motions (e.g. see 
the discussion in |Marconi et al. 2006| . In Fig. [8] we show the 
spectroastrometric curves for models differing only in intrinsic 
velocity dispersion. The comparison with Fig. [7] clearly shows 
that, when adopting the same values of cr , no significant differ- 
ence is found regardless of whether the source of line broadening 
is poor spectral resolution or an intrinsic velocity dispersion in 
the source. 

In general, when spurious line broadening is present, the 
high velocity points in the spectroastrometric curve are not en- 
tirely due to the gravitational potential of the BH, and the curve 
is artificially stretched in the velocity direction leading to a pos- 
sible overestimate of the BH mass with the method described in 
SecE 

From the simulations presented in Fig.[7]and[8]we can verify 
that the effect of line broadening on spectroastrometric curves is 
indeed approximable as a x-axis "stretching". Following equa- 
tion|2]we verified that can then recover the "de-stretched" veloc- 
ities as: 



lobs 



i + 



cr- 



r 2 \ 
_0 

2 



-1/2 



(3) 



thus correcting the spectroastrometric curve. A detailed analysis 
of the effects of spurious line broadening is beyond the scope of 
this paper and will be presented elsewhere but, briefly, we can 
in principle estimate the spurious line broadening cr and that 
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due to unresolved rotation cr by modeling the classical rotation 
curves (e.g. |Marconi et al.|2006 and references therein) and then 
correct the spectroastrometric curves by de-stretching the veloc- 
ity axis as indicated in eq.[3] Clearly, the most accurate approach 
will be that of a combined fit of the classical rotation and spec- 
troastrometric curve, much beyond the simple analysis presented 
here in Sec.|5] 

3.2.7. Flux distribution 

In this section we present simulations based on simple intrin- 
sic flux distribution models, with the aim of showing how sub- 
resolution variations of the flux distribution influence the spec- 
troastrometric curve. 

The basic model adopted for the flux distribution in these 
simulations is an exponential function I(r) = A<?^''° where r is 
the radial distance from the symmetry center located at the po- 
sition (xo,yo) = (0",0") and ro is a characteristic radius which 
we assume equal to 0.05". At the spatial resolution of our simu- 
lations (FWHM =! 0.5"), this emission is spatially unresolved. 




Velocity (km/s) 



Fig. 9. Comparison of spectroastrometric curves for spectra dif- 
fering only in intrinsic flux distribution. Solid line: base model 
flux distribution (explained in the text). Dotted line: model with 
ro = 0.1". Dashed line: model with ro = 0.15". Dot-dashed line: 
base model with a central hole of radius r/, = 0.04". 

In the first set of simulations, shown in Fig. [9] we vary the 
characteristic radius ro in size and introduce a central hole in the 
flux distribution with radius r/,. In the second set of simulations, 



shown in Fig. 10 we vary the position of the center of the flux 
distribution. 

In the "high velocity" range the variation in the centroid po- 
sition with respect to the base model is lower than 0.06" and 
0.09" for the first and second set of simulations, respectively; 
these values should be compared with the pixel size, 0.125", 
and the spatial resolution, FWHM = 0.5". As we will show 
in Sec. [5] these variations will produce only small changes (up 
to +0.2 dex) on the BH mass which can be inferred from spec- 
troastrometric data. On the contrary, much larger variations are 
present in the "low velocity" range. 

This behavior, observed also in previous simulations, is ex- 
plained by the fact that that the gas with lower velocities is lo- 
cated in an extended (spatially resolved) region, and the light 
profiles along the slit are influenced by the shape of the flux dis- 




Velocity (km/s) 



Fig. 10. Comparison of spectroastrometric curves for spectra dif- 
fering only in intrinsic flux distribution. Solid line: base model 
flux distribution (explained in the text). Dotted line: base model 
centered on (0.05", 0.05") position. Dashed line: base model 
centered on (0.1", 0.1") position. Dot-dashed line: base model 
centered on (0.15", 0.15") position. 



tribution. Conversely the gas at higher velocities is confined in a 
small (spatially unresolved) region close to the BH and the light 
profiles along the slit are not influenced by the shape of the flux 
distribution. 

Only in the case of a central hole in the flux distribution does 
the largest variation in the spectroastrometric curve occur in the 
high velocity range. This behavior is explained by the fact that 
that the gas with higher velocities (located in a small region close 
to the BH) is not illuminated because of the central hole of the 
flux distribution and the centroid position is influenced by the 
flux distribution at larger distance from the BH. 

3.2.8. Summary of simulation results 

We can now summarize the results from the above simulations 
and the considerations that can be derived. 

- The presence of a black hole is revealed by a turn-over in 
the spectroastrometric curve, with the high velocity compo- 
nents approaching a null spatial offset from the center of the 
galaxy. Conversely, in the case of no black hole, the curve 
shows a monotic behavior. 

- The information about the BH is predominantly encoded in 
the "high velocity" range of the spectroastrometric curve 
which is generated by spatially unresolved emission. 

- The "high velocity" range of the curve is not influenced 
strongly by the spatial resolution of the observations, leav- 
ing the BH signature unaltered. According to our simulations 
we are able to detect a BH whose apparent size of the sphere 
of influence is as small as ~ 1/10 of the spatial resolution. 
Instead, for the "standard" rotation curves method the infor- 
mation about the BH is effectively canceled when the sphere 
of influence of the BH is not spatially resolved. 

- The amplitude of the spectroastrometric curve decreases by 
increasing the angle between the slit and the line of nodes of 
the gas disk, according to a cosine law. 

- The slit width affects the spectroastrometric curve, due to 
the truncation of the iso-velocity contours. This effect can be 
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reduced to a negligible level selecting for the observations a 
slit width smaller than (or equal to) the spatial resolution. 

- The effect of the finite spectral resolution is to strech (in ve- 
locity space) the spectroastrometric curve. A robust BH mass 
estimates requires a velocity resolution smaller than the line 
width due to (spatially) unresolved rotation. 

- Variations of the intrinsic line flux distribution at sub- 
resolution scales, do not greatly affect the spectroastrometric 
curve in the high velocity range. 

4. Practical application of the method 

4. ) . The use of multiple spectra of the source: the 
spectroastrometric map 

In order to improve the constraints on model parameters with the 
"standard method", it is common practice to obtain several long- 
slit spectra at different slit orientations. Also spectroastrometry 
benefits from this approach as it is possible to recover for the 
same source several spectroastrometric curves, one for each po- 
sition angle (PA) of the slit. With respect to the standard rota- 
tion curves method, where spectra taken along paraller slits have 
been often used, here the fundamental requirement is that the 
various slits must be centered on the expected BH position. 

We simulated the case of three longslit spectra, oriented par- 
allel (PA1), perpendicular (PA2) and forming an angle of 45° 
with the line of nodes (PA3). In Fig. 1 1 we can see the spec- 
troastrometric curves obtained from the three spectra. 



straight line which identifies the line of nodes. Indeed, whenever 
the intrinsic distribution of line emission is circularly symmet- 
ric and centered on the BH position (as assumed in our simu- 
lations), the light centroids at a given velocity are located on 
the disk line of nodes. In the right panel of Fig. 12 we show the 



three-dimensional representation of the spectroastrometric curve 
combining the 2D spatial map of the left panel with the velocity 
axis. In the 3D representation all the points lie on a plane paral- 
lel to the velocity axis and aligned with the disk line of nodes. 
On that plane, the spectroastrometric curve present the usual S 
shape symmetric with respect to the zero point of velocities. 

We remark that in this paper we only apply the spectroas- 
trometric method to continuum subtracted spectra. When deal- 
ing with real data an accurate continuum subtraction will be an 
important requirement to obtain an accurate spectroastrometric 
curve, especially for the "high velocity" points where the line- 
continuum contrast is lower. Any residual continuum emission 
can alter the spatial profile of the emission line, and the "turn 
over" signature of a point mass in the spectroastrometric curve. 
In the latter case, the spectroastrometric displacement of the 
emission line will be shifted towards the continuum photocen- 
tre in the line wings no matter the actual spatial distribution of 
the high velocity gas. The accuracy of continuum subtraction, as 
well as the signal to noise ratio at high velocities, is therefore 
a fundamentally limit to the ability of detecting the central BH. 
A more detailed analysis of the effects of continuum emission 
will be presented in a forthcoming paper where we will apply 
the spectroastrometric method to real data. 
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Fig. 11. Spectroastrometric curves for three simulated spectra of 
the same model with different PA of the slit. Solid line: slit with 
PA (referred to the disk line of nodes position angle) 6 s iit~SwN = 
0°. Dotted line: slit with e sUt - LO n = 45°. Dashed line: slit with 
Qsttt - Own = 90°. 

Each spectroastrometric curve provides the photocenter po- 
sition along one slit, i.e. the position of the photocenter projected 
along the axis identified by the slit. Combing the spectroastro- 
metric curves we can thus obtain the map of photocenter posi- 
tions on the plane of the sky for each velocity bin. In principle, 
the spectroastrometric curves from two non-parallel slits should 
suffice but we can use the redundant information from the three 
slits to recover the 2D sky map as described in more detail in 
appendix [B] 

In Fig. 12 (left panel) we present the 2D spectroastrometric 
map on the plane of the sky. All the points of the map lies on a 



4.2. The effect of noise on spectroastrometry 

We now consider the effect of noise on spectroastrometric curves 
and maps, as in real data. We simulate noise by adding normally 
distributed random values to the simulated position-velocity 
spectral maps. These random numbers are characterized by zero 
mean and standard deviation cr noise which is chosen in the fol- 
lowing way. We first extract a "nuclear" spectrum by coadding 
line emission over an aperture equal to the PSF FWHM, cen- 
tered on the nucleus position. The S /N ratio of the spectrum is 
then defined on the peak of the line profile, and the actual cr no i se 
value per pixel is selected to provide a given S /TV. 

In Fig. 13 we show the spectroastrometric curves and maps 
for the same model presented in Fig. TT|to which we added an 



artificial noise in order to obtain S/N of 20, 50 and 100, respec- 
tively. 

A key point to emphasize is that the curves (and maps) with 
higher S/N tend to have more measured points. This is because at 
the "high velocities" we have a low flux with respect to the line 
core and the presence of noise limits our ability to obtain a reli- 
able light centroid position. In detail we stop estimating the light 
centroids when the mean flux of the light profile for that partic- 
ular velocity bin is lower than the noise level. Thus, a spectrum 
with higher S/N enables us to obtain more accurate measure- 
ments on individual data-points, but also extends the coverage at 
higher velocities, crucial to detect a BH and to measure its mass. 

The 2D spectroastrometric maps can be used to estimate the 
geometrical parameters of the gas disk. In the following, we will 
always make use of simulations with noise. 
As noted in Sect. 



3.2.1 



if the gas kinematics is dominated 
by rotation around a pointlike mass (the BH), the position of the 
light centroid at the highest velocities approaches the position 
of the BH. This consideration, which is also valid for the 2D 
spectroastrometric map, allows us to estimate the position on 
the plane of the sky of the BH as the average position of the 



Gnerucci et aL: Spectroastrometry of rotating gas disks and supermassive black holes. 



11 




Fig. 12. Spectroastrometric 2D map derived from the three spectroastrometric curves of Fig. [nj Left panel: derived photocenter 
positions on the sky plane, the black open diamonds are the points actually used for the y 2 minimization, the dashed lines indicate 
the three slits each with its derived centre location (filled squares). Right panel: the 3D plot of the map, where the z axis is velocity. 



Table 1. Lines of Nodes and BH position estimates from noisy 
data 



Model 


x BH (") 


)>bh(") 


OlonC) 


Noiseless 


0.00 


0.00 


0.00 


S/N=20 


-0.01±0.02 


0.00±0.03 


-5±4 


S/N=50 


0.00±0.01 


0.00±0.01 


-1±1 


S/N=100 


0.015±0.008 


-0.003±0.008 


-0.9±0.8 



"high velocity" points. In practice, we calculate this position by 
first taking the averages of the coordinates of the points in the 
blue and red "high velocity" ranges and then taking the average 
coordinates of the "blue" and "red" positions. The inferred BH 
position is marked by a red filled circle in the right row diagrams 
of Fig. [13 



Also, the high velocity points in Fig.[l3](and also Fig.[L2} lie 
on a straight line which identifies the direction of the disk line 
of nodes (see Fig.|6]and Sect. 3.2.4 1. The red solid line shown 
in the right row of Fig. |4] represents the position angle of the 
line of nodes (9[£>n) obtained by fitting a straight line to the high 
velocity points. 

In the Table [T] we report the BH position and disk line of 
nodes position angle values obtained from the 2D map of the 
model of Fig. [13] We can observe that the accuracy of these mea- 
surements scales approximately with the square root of the S/N 
of the spectrum. Moreover in the maps with lower S/N the scat- 
ter of the points around the line of nodes increases (the points 
are almost perfectly aligned on the line of nodes direction in the 
noise free model of Fig. 10 1. This fact decreases the accuracy of 
the BH position and line of nodes PA estimates for lower S/N. 



5. Estimate of the BH mass from the 
spectroastrometric map 

In this section we present a simple and straightforward method 
to recover the BH mass from the spectroastrometric measure- 
ments. The method is based on the following assumptions: (i) 
the gas is circularly rotating in a thin disk and (ii) each point of 
the spectroastrometric map represents a test particle on the disk 
line of nodes and is characterized by a line of sight ("channel") 
velocity which is given by the center of the velocity bin where 
the light centroid was estimated. Under these assumptions it is 
trivial to relate the spectroastrometric map to the BH mass. 

The circular velocity of a gas particle with distance r from 
the BH is given by: 



Vr, 



G[M BH + M star (r)] 



(4) 



where M slar (r) is the stellar mass enclosed in a sphere of radius 
r (assuming spherical symmetry) and can be written as: 



M star {r) = M/L ■ L{r) 



(5) 



with L(r) representing the radial luminosity density distribution 
in the galactic nucleus and M/L is the mass to light ratio of the 
stars. 

The line of sight velocity (hereafter V c i, for "channel veloc- 
ity") of a test particle located on the disk line of nodes at distance 
r from the BH is then given by 



V c h = V ro ,sin(i) + V s . 



(6) 

is the sys- 



where i is the inclination of the disk plane and V y 
temic velocity of the galaxy. 

The centroid positions on the spectroastrometric map are de- 
noted by (x c h,y c h) and, in general, they are not perfectly aligned 
along the line of nodes as described in Sect. 4. 1 and B] We select 
the "high velocity" points and we then identify the ine of nodes 
by a linear fit of the (x c ;,,y c /,) spectroastrometric map. To reduce 
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Fig. 13. Left row: spectroastrometric curves for the model of Fig. [IT] with added noise for a S/N of 20, 50, and 100 (from top to 
bottom). Solid line: slit at PA = 0°. Dotted line: slit at PA 45°. Dashed line: slit at PA = 90°. For simplicity we plot the error bars of 
the points only for the PA= 0° curve. Right row: spectroastrometric 2D map derived from the three spectroastrometric curves of the 
left row. The red filled circle marks the derived position of the BH, while the red solid line represents the position angle of the line 
of nodes (Own) obtained by fitting a straight line to the X-Y position of the high velocity points. 



the correlation between the values of intercept and slope of the 
fitted line we consider 

y = ( x - x mean ) tan{6 L0N ) + b (7) 



where x mean is the mean of the x C h positions. Olon and b denote 
respectively the position angle and the intercept of the line of 
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nodes to be determined by the fit. Then Qlon and b are recovered 
by minimizing the following x 1 where we take into account both 
errors in x and y: 



X 



-z 

ch 



bch - (x c h - x mean )tan(0 LON ) - b] 
Ay] h + tan(0 WN ) 2 Ax 2 ch 



(8) 



where Ax c i wn and Ay c i, an are the uncertainties on the position 
of the points of the spectroastrometric map and, as observed in 
Sect. 4.1 the sum is extended only over the "high velocities" 
range. 

The second step requires to associate each velocity bin to its 



position within the rotating disk. As described in Sect. 4. 1 the 



2D map points lie with good approximation on the line of nodes, 
particularly in the "high velocities" range. We then project the 
(x c h,y c h) position of the 2D-map points on the line of nodes and 
calculate their coordinate "S" with respect to this axis of refer- 
ence, with "S" defined as: 



Sch= X ch COs(0 LO N) + \y, 



ch 



tan{6 LO N)]sin{0 LON ) (9) 



Each S c h has an uncertainty AS c y, estimated by taking into ac- 
count the uncertainties on the 2D map points position (Ax C f, and 
Ay c h as well as the uncertainties on the line of nodes parameters 
(A6lon, Afe) resulting from the fit. The distance r of a test parti- 
cle from the BH used in Eq.|4]differ from S c h only by a constant 
So, the unknown coordinate of the BH along the line of nodes, 
i.e. 



r c h = \S 



ch 



(10) 



Combining Eqs. |4J [5] |6J and [TO] we finally obtain the expected 
velocity in channel ch: 



Vrh = 



G[M BH + M/L-L(\S ch -S \)] 



\S ci, - S c 



sin(i) + W 



(11) 



The unknown model parameters are: 



M BH mass of the BH 

Ml L mass to light ratio of the nuclear stars 
So line of nodes coordinate of the BH 
Vsy S systemic velocity of the galaxy 
i inclination of the gas disk 

and can be determined by minimizing the following x 1 

2 



X 



-z 



Wch - V sys \ - \V ch {par) - V sys \ 



A(S ch ;par) 



(12) 



where A(5 e /,; par) is the uncertainty of the numerator com- 
puted as a function of the unknown parameters values (par) and 
the uncertainties AS c h on the positions along the line of nodes. 
We remark that the channel velocity V c h has no associated uncer- 
atinty since it is not a measured quantity but is the central value 
of a velocity bin where the spectroastrometric curve is measured. 
As discussed in the previous sections, we restrict the fit (i.e. the 
sum over the velocity channels) to the "high velocities" range. 

is much smaller for 



The A(S c h\ par) factor in equation 12 



the points at lower velocities (i.e. smaller \V c h - V sys \) which are 
closer to the peak of the line profile and have therefore much 
larger S /N than the points at "high velocity". However, from the 
discussion in the previous sections, we know that the spectroas- 
trometry is less reliable for the determination of BH properties 
like position and mass. Therefore, to avoid being biased by po- 
tentially faulty points we add in quadrature a constant error A sys 
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Fig. 14. Fit results for models of Fig. 
panel), S /N = 50 (middle) and S IN 
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with S/N = 20 (top 



100 (bottom). Points with 
error bars are the result of "observations". The red line is the 
best fit model and the dotted line is the expected contribution 
from the stellar mass. 



to A(5 C /,; par). The value of A sys is found by imposing that the 
reduced y 1 is equal to 1 . Therefore, if the points at lower veloc- 
ity with higher S /N are problematic they will not bias the final 
fit results, because their high weight will be greatly reduced by 
the addition od of a much larger A sys . 
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Fig. 16. Effect of varying the intrinsic flux distribution at sub-resolution scales: results for various intrinsic flux distributions. Top 
left panel: exponential flux distribution centered on (0", 0") with ro = 0.05". Top rigth panel: exponential flux distribution centered 
on (0", 0") with ro = 0.15". Bottom left panel: exponential flux distribution centered on (0", 0") with ro = 0.05" and a central hole 
of radius r/, = 0.04". Bottom rigth panel: exponential flux distribution centered on (0.1", 0.1") with ro = 0.2". 



After determining the best fit values of the free parameters, 
the position of the BH on the plane of the sky is simply given by 



Xbh = Sqcos{9 LO n) 

ytiH = S sin(6 L0 N) + b- x„ 



1 

tan(6 L0N ) 



(13) 
(14) 



We have performed the fit of the spectroastrometric data for the 
models of Figs. 13 which are computed with the same param- 
eters but different spectrum S IN. The results of the fitting pro- 
cedure are presented in Fig. 14 where we plot the high veloc- 
ity points of the spectroastrometric map in the \V c f, - V sys \ vs 
r (|S c h - Sol) diagram. The solid red lines represent the curves 
expected from the rotating disk model. In Fig. 
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we present 

the direct comparison between the observed spectroastrometric 
curve along the line of nodes (e.g. S c h im vs V c han) with the best 
fit model (for S /N = 100). 

In all cases, for each "test particle", we are able to measure 
photocenter positions on the plane of the sky with an accuracy 



of better than 07025, corresponding to 1/20 of the spatial res- 
olution (FWHM = 075). 

The comparison between best fit and simulation parameters 
is shown in Table [2] In particular, we can recover the BH mass 
value with an accuracy of ~ 0.1 dex. The M/L value is not well 
constrained by the models because with the adopted parameters 
the contribution to the gravitational potential of the stellar mass 
is negligible. 

In the test cases considered here (Mbh = 10 Mq, stellar 
velocity dispersion of cr„ flr - 200 km/s and D 3.5Mpc), the 
radius of the BH sphere of influence is r BH =s 10 pc correspond- 
ing to =s 076, barely resolved with the adopted 075 resolution. 
With the spectroastrometric technique, we have instead probed 
down to radii of the order of 0.05" (~ 1 pc), a factor ~ 10 smaller 
than rBH- This clearly demonstrates that with spectroastrometry 
we can recover information on spatial scales which are much 
smaller than the limit imposed by the spatial resolution, thus 
opening the possibility of probing the gravitational potential in- 
side smaller BH spheres of influence. 
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Fig. 17. Application of the method to low BH masses: results from simulations with S/N = 100, M BH = 1O 7 M0 (left) and 10 6 5 M© 
(right) . Top: spectroastrometric curves. Middle: 2D spectroastrometric maps. Bottom: fit results. 



We now present a set of simulations where we vary the in- 
trinsic flux distribution of the gas at sub-resolution scales and 
observe how this affects the recovery of the input model param- 
eter values from the fit. As in section 3.2.7 we adopt a basic 



model for the flux distribution that is an exponential function 
(I(r) = Ae r l r °) where the symmetry is center located at the po- 
sition (xo,yo) = (0",0") and the characteristic radius is set to 
ro = 0.05". In Figs. 16 we show the cases of: the basic flux dis- 



tribution, flux distribution with characteristic radius ro = 0.15", 
flux distribution with central hole of radius = 0.04", flux dis- 
tribution with characteristic radius ro = 0.05" but centered at the 
position (xo,yo) = (0.1", 0.1"), all with a simulated noise for a 
S/N of 100. In Table [3] we report the best fit values of the free 
parameters. These simulations are chosen to represent extreme 
cases of sub resolution variation of the intrinsic flux distribution. 
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Table 2. Fit results from the baseline model with different S /N. 



Parameter 




Fit result 




Model value 




a/JN = 1 UU 


S/N=50 


S/N=20 




e L0 N [°] (WNfit) 


-0.9 ± 0.8 


-1 ± 1 


-5 ±4 


0.00 


b ["] (WNfit) 


-0.001 +0.002 


-0.002 ± 0.003 


0.000 ± 0.009 


0.00 


log l0 (M BH /M Q ) 


8.05 ± 0.03 


8.00 ± 0.06 


7.9 ±0.1 


8.00 


log w (M/L) [Mq/Lq\ 


— y.Ul 




U.2 ± 2.6 


U.UU 


So ["] 


0.02 ± 0.02 


-0.03 ± 0.03 


-0.02 ± 0.03 


0.00 


Vsy S [km/s] 


528 ± 26 


478 ± 33 


489 ± 23 


500 


i [°] 


35.0 * 


35.0* 


35.0 * 


35.0 " 


A. [VJ [km/s] 


25 


28 







reduced^ 2 (x 2 /D.O.F.) 


1.003 (5.12/5) 


0.992 (2.976/3) 


0.015 (0.031/3) 




*BH ["] 


0.02 ± 0.02 


-0.03 ± 0.03 


-0.02 ± 0.03 


0.00 


JBH ["] 


-0.001 ±0.002 


-0.001 ± 0.004 


-0.002 ± 0.009 


0.0088 



" Parameter not constrained from the fit 
* Parameter hold fixed 



Table 3. Effect of varying the intrinsic flux distribution at sub-resolution scales: fit results. 



Parameter 




Fit result 




Model value 




flux distribution A-' 


flux distribution 


flux distribution C' 


flux distribution D f 




9lon [°] (WNfit) 


-1.0 ± 1.3 


0.9 ± 1.1 


-4.1 ± 1.3 


5.7 ±8.3 


0.00 


b ["] (WNfit) 


0.001 ± 0.002 


0.000 ± 0.002 


0.002 ± 0.003 


0.02 ± 0.01 


0.00 


log m (M BH /M Q ) 


7.91 ±0.02 


8.16 ±0.03 


8.04 ± 0.03 


8.13 ±0.10 


8.00 


log m (MIL) [Mq/Lq] 


-8.23 a 


-7.41 " 


-8.15° 


-6.68 " 


0.00 


So ["] 


0.012 ± 0.008 


-0.01 ±0.03 


-0.01 ± 0.01 


0.02 ± 0.03 


0.00 


V svs [km/s] 


534 ± 18 


494 ± 32 


478 ± 35 


500* 


500 


i [°] 


35.0* 


35.0* 


35.0* 


35.0* 


35.0 


A JYI [km/ s] 


8 


28 


22 







reduced^ 2 (x 2 /D.O.F.) 


1.05 (5.25/5) 


1.05 (4.19/4) 


1.05 (6.28/6) 


1.006 (1.006/1) 




*BH ["] 


0.012 ±0.008 


-0.01 ± 0.03 


-0.01 ±0.01 


0.02 ± 0.03 


0.00 


JBH ["] 


0.001 ± 0.002 


0.001 ± 0.003 


0.001 ± 0.002 


0.02 ± 0.02 


0.00 



f Flux distribution. A: exponential centered on (0",0") with r = 0.05". B: exponential centered on (0",0") with r = 0.15". C: exponential 

centered on (0", 0") with ro = 0.05" and a central hole of radius r;, = 0.04". D: exponential centered on (0.1", 0.1") with ro = 0.05". 
" Parameter not constrained from the fit 
* Parameter hold fixed 



We can observe that with the fit we can recover the correct 
BH mass value with an accuracy of better than ~ Q.\6dex, show- 
ing that even "extreme" variations of the flux distribution at sub- 
resolution scales, do not greatly affect the BH mass estimate. 

Finally we consider a set of simulations where we decrease 
the M BH value. In Figs. 17 we show the cases of BH masses 



M B h = 10 7 Mq and M BH = 10 6 - 5 M© respectively, both with a 
simulated noise for a S/N of 100. In Table [4] we report the best 
fit values of the free parameters. 

The case of lower black hole mass is particularly significa- 
tive. We choose a M BH = 10 6 5 Mq value because, as observed in 
section [3. 2. 8 [ the standard rotation curve method cannot distin- 
guish between a case of 1O 6 5 M0 BH and no BH. By looking at 
Fig. 17 we note that the BH signature is still present in a realistic 



spectroastrometnc curve. 

Nonetheless, the fit is not well constrained by the few dat- 
apoints. Furthermore, we obtain a value of M B h of \Q 12 Mq 
(Table substantially larger than the model mass of 10 6 ' 5 M o . 



This is due to the fact that in our simulation the size of the 
velocity bins in the PVD diagram is 50 km/s, matched to the 
assumed spectral resolution (R ~ 6000). This value is insufficient 
to adequately resolve the emission line profile. As observed in 
sec. 3.2.5 the spectroastrometric curves are stretched along the 



velocity axis with artificially high velocity values which result 
in an overestimated BH mass. 

This problem can be overcome by increasing the spectral res- 
olution of the simulated data to a higher value, R ~ 30000, still 
well within reach of existing spectrographs. We then repeated 
the simulation setting the size of the velocity bin of the longslit 
spectra to 10 km/s. In Fig. [18] we show the spectroastrometric 
curves and the result of the fit for this model. The sampling of 
the spectroastrometric curves is largely improved and the fit is 
consequentely better constrained. The derived BH mass value, 
\og(M BH / Mq) = 6.8 ± 0.1 is now in better agreement with the 
model value. The results of these simulations confirm that to 
assess the accuracy of BH masses from the spectroastrometric 
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Parameter 


M B h = 


10 7 M o model 




M BH = 10 6 5 M o model 






Fit result 


Model parameters 


Fit result 


Fit result B Model parameters 


Olon [°] {line of nodes fit) 


+ 2 


0.00 


-2 ±5 


-0.5 ± 0.7 


0.00 


b ["] (line of nodes fit) 


-0.007 ± 0.005 


0.00 


-0.02 ± 0.02 


-0.001 + 0.002 


0.00 


log w (M BH ) [log w (M Q J] 


7.20 ± 0.03 


7.00 


7.23 ± 0.04 


6.82 ± 0.08 


6.50 


log l0 (M/L) [log l0 (MQ/LQ)\ 


n nn b 
U.UU 


n on 
U.UU 


o on b 
U.UU 


0.00* 


0.00 


So ["] 


-0.03 + 0.02 


0.00 


0.00 * 


0.025 ± 0.00 c 


0.00 


V m [km/s] 


493 ±6 


500 


500* 


503 ±4 


500.0 


i [°] 


35.0 6 


35.0 


35.0* 


35.0* 


35.0 


Ajyj [km/s] 










8 




reduced chi 2 (x'/D.O.F.) 


0.04 (0.04/1) 




0.26 (0.26/1) 


1.13 (1.13/1) 




XBH ["] 


-0.03 ± 0.02 


0.00 


0.00 ± 0.00 * 


0.025 ± 0.00 d 


0.00 


}'BH ["] 


-0.007 + 0.005 


0.00 


-0.02 ± 0.02 


-0.001 ± 0.003 


0.00 



A Simulation of the model with 50km/ s velocity bins 

B Simulation of the model with 10km/ s velocity bins 

* Parameter hold fixed 

c Parameter at the edge of allowed range 

d Parameter calculated from fixed parameters 
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Fig. 15. Observed and model spectroastrometric curves at S /N — 
100 (e.g. Schan vs V c i mn ). The horizontal dotted line is the BH "S" 
coordinate So- The vertical dotted line is the systemic velocity 



method one should first check that the line widths in the central 
emission region are well resolved spectrally. 

In the last fits reported in table |4] we hold fixed the mass- 
to-ligth ratio M/L to the model value. The number of "high ve- 
locity" points in the spectroastrometric curves is small (only 2 
for the model M BH = 1O 6 5 M0 with 50km/s velocity bin) and 
it is not possible to constraint the radial dependence of the mass 
distribution. In general, when dealing with real data with very 
few "high velocity" points, one can only measure the total mass 
enclosed in the smaller radius which can be estimated from spec- 
troastrometry. If it is possible to analyze the "classical" rotation 
curve, one could determine the mass-to-ligth ratio M/L from that 
and use it in the spectroastrometric analysis. The fixed mass-to- 
ligth ratios used in the above analysis should be interpreted as re- 
sulting from a classical analysis of the extended rotation curves. 



6. Integral Field Spectroscopy 

The extension of the spectroastrometric technique from multi- 
slit observations to Integral Field Units (IFU) is straightforward, 
and carries with it many advantages. With IFUs the complex is- 
sue of the limited spatial coverage of the slits (leading to the 
problem of velocity truncation) is partially removed and a full 
2D map can be directly derived. Another substantial advantage 
is the reduced observing time requested to obtain a given S/N 
level, since there is no need to obtain observations of the same 
galaxy at different silt position angle. 

The analysis of IFU data now reduces to fitting a 2D gaus- 
sian to each channel map in turn yielding the X,Y photocenters 
as a function of velocity, i.e. the spectroastrometric map. Armed 
with these derived photocenters one proceeds directly to the ap- 
plication of Sect. [5] 

In practice it is important to determine the accuracy of each 
photocenter. This can easily be accomplished by using monte 
carlo realizations drawn from the original data in each channel 
giving a distribution of i values for X and Y from which one can 
determine a median and inter-quartile spread. 

In order to provide a quick analysis of advantages and disad- 
vantages in the use of IFUs ve longslit spectra, we first recall the 
mandatory requirements needed to obtain an accurate and useful 
spectroastrometric curve. The spectra must be characterized by 
good spectral resolution, signal-to-noise ratio and spatial sam- 
pling. The emission line profile must be well spectrally resolved 
otherwise, as discussed in Sec. |3. 2. 5\ one would overestimate the 
BH mass. High signal-to-noise and spatial sampling are required 
for a robust and accurate robust estimate of the spatial centroid 
of the emission line, as discussed in detail in appendix [A] 

Compared to longslit spectrographs, IFUs have the advan- 
tage of a two-dimensional spatial covering of the source which 
provides a more direct and accurate determination of the spec- 
troastrometric map on the plane of the sky. However, for a given 
number of detector pixels, this is usually done at the expense 
of spectral resolution thus limiting the application of the spec- 
troastrometric analysis to larger BHs. As a compromise between 
field-of-view and spectral resolution, an IFU can have a lower 
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Fig. 18. Application of the method to low BH masses: results 
from the simulations with BH mass of Mbh = 10 6 - 5 M© and ve- 
locity bins of lOkm/s. Top: Spectroastrometric curves (for sim- 
plicity we plot the error bars of the points only for the PA= 0° 
curve). Middle: spectroastrometric map. Lower panel: The result 
of the fit. 



spatial sampling compared to longslit spectra which limits the 
accuracy of the photocenter determination. However, with IFUs 
it is possible to integrate longer, since one does not need to ob- 
tain spectra at different position angles. 



As an example, let us compare the seeing-limited perfor- 
mances of ISAAC and SINFONI (both at the ESO VLT) for 
spectroastrometric use. ISAAC has a spatial sampling along the 
slit of 0.147", where for SINFONI in the seeing limited mode 
spatial pixels (spaxels) have angular dimensions of 0.125" x 
0.25". There exist modes with finer spatial sampling (0.05" x 
0.1") but at the expense of a smaller field of view which might 
be insufficient in seeing limited observations. Overall, the spatial 
sampling in ISAAC and SINFONI are not extremely different, 
at least for moderately good seeing (~ 075). Regarding the spec- 
tral resolution, let consider the case of the K (J) band, at 2.2fan 
(~ 1 2fxm). ISAAC, in the Short Wavelength Medium Resolution 
configuration offer a spectral resolution of 8900 (10500) while 
SINFONI with the K (J) grating offers a spectral resolution of 
4000 (2000). Clearly, the advantage of the 2D coverage provided 
by SINFONI is obtained at the expense of spectral resolution, 
and therefore SINFONI is not a good choice to detect small BHs 
(or large BHs in more distat galaxies). 

In conclusions, IFUs have the advantage of the two- 
dimensional spatial covering which allows longer integrations 
on source and more accurate spectroastrometric maps. On the 
other hand, longslit spectrographs can provide a better combi- 
nation of high spectral resolution and spatial sampling for the 
detection of smaller BHs. 

Regardless of the use of IFUs or longslit spectrographs, the 
super-resolution provided by spectroastrometry allows BH mass 
determinations in galaxies at distances substantially higher than 
those that can be studied with the "standard" rotation curves 
method. 



7. Combining Spectroastrometric and Rotation 
Curves 

In this paper we have shown that by means of spectroastromet- 
ric curves we can recover information at scales smaller than the 
spatial resolution of the observations and we have provided a 
simple method to estimate BH masses. However, spectroastrom- 
etry is not a replacement of the standard rotation curve method. 
As shown in Sec. 3 (e.g. Fig. 1), the spectroastrometric curve is 
complementary to the rotation curve and it is clear that any kine- 
matical model must account for both curves at the same time. 
For example, any contribution from extended mass distributions 
(e.g. stars) can be well constrained with rotation curves but not 
with spectroastrometric curves which sample very small scales. 
A detailed analysis of the constraints posed simultaneously by 
spectroastrometric and rotation curves will be discussed in forth- 
coming papers. 

The reader might wonder why a combined modeling of spec- 
troastrometric and rotation curves should be better than than 
modeling the full Position- Velocity diagram. The reason is very 
simple and is related to the fact that the full PVD depends on the 
unknown intrinsic flux distribution of the emission line. M arconil 
|et al.| ( [2006| showed that line profiles do depend on the assumed 
flux distribution while, on the contrary, such dependence is much 
weaker on the first moment of the line profile (i.e. the mean ve- 
locity). In this paper we have shown that the spectroastrometric 
curve depends on the assumed flux distribution in the low veloc- 
ity range, while such dependence almost disappears in the high 
velocity range. Therefore, by selecting the rotation curve and 
the spectroastrometric curve we can greatly diminish the effects 
of the unknown intrinsic flux distribution which plague the full 
PVD. 
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8. Summary and Conclusions 

In this paper we have discussed the application of the spectroas- 
trometric method in the context of kinematical studies aimed at 
measuring the masses of supermassive black holes in galactic 
nuclei and we have presented its advantages compared to classi- 
cal "rotation curves" method. 

We have conducted an extensive set of simulations which 
have shown that the presence of a supermassive black hole is re- 
vealed by a turn-over in the spectroastrometric curve, with the 
high velocity component approaching a null spatial offset from 
the location of the galaxy nucleus. All the relevant information 
about the BH is encoded in the "high velocity" range of the spec- 
troastrometric curve, which is almost independent of the spatial 
resolution of the observations. According to our simulations, the 
use of spectroastrometry can allow the detection of BHs whose 
apparent size of the sphere of influence is as small as ~ 1 / 10 of 
the spatial resolution. 

We have then provided a simple method to estimate BH 
masses from spectroastrometric curves. This method consist in 
the determination of the spectroastrometric map, that is the po- 
sitions of emission line photocenters in the available velocity 
bins. This is trivially obtained from IFUs but can also be ob- 
tained by combining longslit spectra centered on a galaxy nu- 
cleus and obtained at different position angles. From the "high 
velocity" points in the spectroastrometric map one can then ob- 
tain a rotation curve to trivially estimate the BH mass. We have 
provided practical applications of this method to simulated data 
but considering noise at different levels. From this analysis, we 
confirm that with seeing limited observations (~ 075) we are 
able to detect a BH with mass 10 6 5 M Q (D/3.5 Mpc), where D 
is the galaxy distance. This is a factor ~ 10 better than can be 
done with the classical method based on rotation curves. 

Finally, we have discussed advantages and disadvantages of 
IFU vs longslit spectrographs concluding that IFUs have the ad- 
vantage of the two-dimensional spatial covering which allows 
longer integrations on source and more accurate spectroastro- 
metric maps. On the other hand, longslit spectrographs can pro- 
vide a better combination of high spectral resolution and spa- 
tial sampling for the detection of smaller BHs. Regardless of the 
adopted type of spectrograph, the super-resolution provided by 
spectroastrometry allows BH mass determinations in galaxies at 
distances substantially higher than those that can be studied with 
the "standard" rotation curves method. 
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Appendix A: The measure of the light profile 
centroid 

The first fundamental step in the application of the spectroastro- 
metric method is to find the best possible way to estimate the 
centroid position of the light profile for a given wavelength (or 
velocity). In principle one could simply estimate the centroid po- 
sition as the weighted mean of the position with the correspond- 
ing emission line flux. If y, denotes the position and /, the flux 
of the i-th pixel along the slit (for a given velocity), the centroid 
position y cent is simply: 



ycent — 



(A.l) 



The problem in using equation |A.1| resides in the presence of 
noise and in the fact that the value of y cenl depends on the choice 
of the position range used to compute the weighted mean. In Fig. 
A. 1 we exemplify the effect of choosing different ranges to cal- 
culate the centroid of real data. This effect will not be negligible 
especially because one is looking for precisions which are much 
lower than the spatial resolution (e.g., ~ 2.5 pixels for the data 
shown in the figure). 



Position profile at v=410km/s 




160 170 180 190 

Position along the slit (pixel) 

Fig. A.l. Example of a light profile along the slit for a given 
velocity bin for real data. Three different position ranges are in- 
dicated as horizontal lines and the corresponding centroid posi- 
tions are indicated with vertical lines with the same style. 

Another possibility to estimate the centroid position is pro- 
vided by a parametric fitting of the light profile with an appro- 
priate function. For instance, we have tried a Gauss-Hermite ex- 
pansion: 



F(x) = Ae T [1 + h 3 H 3 (y) + h 4 H 4 (y)] 



(A.2) 



where y = (x - xq)I<t and H 3 and Ha, are the Gauss-Hermite 
polynomials of the third and fourth order. 

The advantage of the fitting method is that with the paramet- 
ric fit the results are much less sensitive to noise and choice of 
the position range for the analysis but it is difficult to find an 
appropriate parametric function which can reproduce the profile 
shape without a large increase of the free parameters which, in 
turn, decreases the accuracy of the centroid position. The light 
profile along the slit has often a complex shape, sometimes with 
multiple peaks and there is not any physical reason to suggest a 
particular parametric function. 
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Moreover, the shape of the light profile for a given veloc- 
ity bin is determined not only by velocity field of the rotating 
gas disk but also by many other factors, like the intrinsic light 
distribution of the emission line, the slit width and position and 
many others. Therefore one has to investigate what is the signa- 
ture which is most related to the velocity field. 

The piece of information we seek is the average position of 
the gas rotating at a given observed velocity which can be con- 
verted into the gravitational potential in the nuclear region, with 
the usual assumption of a thin, circularly rotating disk. 

As explained above, the use of the entire light profile along 
the slit for a given velocity bin to calculate the photocenter po- 
sition is extremely complex, moreover such position might be 
influenced by many features which do are not connected with 
the rotation of the disk. After testing the method with real and 
simulated data, we have reached the conclusion that the quantity 
we need to measure is not the centroid of the whole light profile 
but only of the central spatially unresolved component which is 
related to gas located very close to the galaxy nucleus (e.g. gas 
rotating around the BH). Thus, the problem of measuring the 
light centroid along the slit can be transformed into measuring 
the position of the peak of the central unresolved emission. 
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Fig. A.2. Example of centroid determination on a real spectrum. 
The dashed line denotes the fitted gauss-hermite functions and 
the relative centroid position (vertical dashed line). The solid 
lines denotes the gaussian function which has been fitted in a 
window around the principal peak as described in the text. The 
vertical solid lines represent the edges of the window and the 
gaussian center respectively. 



The profile peak at the location of the nucleus is due to 
the spatially unresolved emission from the innermost parts of 
the rotating gas disk. Because it is spatially unresolved, the 
width of central peak is expected to be of the order of the spa- 
tial resolution, and its shape can be well approximated with a 
gaussian function. This component is easy to identify and fit. 
Moreover this unresolved component is emitted from the inner- 
most parts of the gas disk, where the gravitational potential of 
the BH is strongest with a larger influence on the gas kinemat- 
ics. Conversely, the spatially resolved features are very complex 
to fit and are influenced by many factors not directly connected 
with the gravitational potential. 

In practice, we first obtain an initial guess of the principal 
peak position and we select the light profile within a box cen- 
tered on this initial guess and whose width is of the order of the 
spatial resolution. We then accurately locate the peak position 
by fitting a simple gaussian function to this selected portion of 
the light profile. In details, we usually set the width of the box 
as ko-Q where cr is the spatial resolution (standard deviation of 
the gaussian approximating the point spread function) and k a 
constant between 1 and =s 3. Using values of k too close to 1 
results in too few fit points but k should also not be too large if 
we want to consider only the unresolved component. From our 
tests, the best choice has come out to be k - 2. 

This method used to obtain the light centroid provides us 
with an important control parameter, that is the width of the fitted 
gaussian: if this is much larger than the spatial resolution, this 
is an indication that line emission is dominated by a spatially 
extended component and the obtained value of the centroid is 
less reliable. 
In Fig. 



A.2 



we show an example of centroid determination 
obtained by fitting a gauss-hermite function to the whole profile 
and a simple gaussian in a 2<x wide window centered around 
the principal peak. The centroid values obtained with the two 
methods might or might not be consistent but it is clear that the 
centroid value with the first method is shifted to the right by the 
contribution of the profile wings, which, moreover, are not even 
well reproduced by the gauss-hermite function. In the bottom 
panel of Fig. A.2 none of the two methods can well reproduce 
the profile peak, but the width of the fitted gaussian turns out 
to be ~ 2.5 times larger than the spatial resolution indicating 
that the principal peak is partially resolved, and the value of the 
centroid is less reliable. 

Following these considerations, we found that the most ro- 
bust and objective way of defining the "high velocity" range is 
through the width of the gaussian fitted to the principal peak 



of the light profile: in Fig. A. 3 we display the FWHMs (Full 



Width Half Maxima) of the gaussian fitted to the principal peak 
of the light profile for spectroastrometric curves of Fig. [TT] We 
can observe that at low velocities the FWHM increases because 
the peak is no longer unresolved. We compare the FWHM to the 
spatial resolution (~ 0.5") because the FWHM of the light pro- 
file of an unresolved source should be of the order of the spatial 
resolution. We selected the "high velocity" range by imposing 
that the FWHM is lower than 1.1 times the spatial resolution 
(FWHM of the PSF) resulting in the range v < -250 km/s and 
v > 250 km/s for the spectroastrometric curves. 

Appendix B: Building the spectroastromertic map 
from multiple slits observations. 

In this section we describe how to build a 2D spectroastrometric 
map, i.e. as to estimate the photocenter position on the plane of 
the sky for each velocity bin, from multiple slits observations. 
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Fig.A.3. FWHM of the gaussian fitted to the principal peak 
of the light profile principal peak for spectroastrometric curves 
of fig. 11 Solid line: PA1 curve. Dotted line: PA2 curve. Dot- 
dashed lines: PA3 curve. The horizontal dashed line denotes 1.1 
times the PSF FWHM. The gray dashed lines are the estimated 
limits of the "high velocities". 



N — 4, where 4 is the number of unknown slit center coordinates 
and 2N is the number of unknown photocenter positions. Since 
N is larger than 4 the problem is well posed. The final spectroas- 
trometric map on the plane of the sky is that given by the best 
fitting set of slit centers. 



The sum over v, in ( |B.3[ ) is extended only over the "high ve 
locity 



range. Indeed, in Sect. |3.2| we concluded that the "high 
velocity" range of the spectroastrometric curve is more robust, 
and less affected by slit losses which artificially change the pho- 
tocenter position in the low velocity range. As explained in ap- 
pendix [A] where we discuss in detail, centroid determination, 
we find the spectroastrometric curve by fitting a gaussian to the 
principal peak of line emission along the slit (that centered on 
the nucleus position). Then the "high velocity range" is selected 
by imposing that the FWHM of the fitted gaussian is lower than 
1.1 times the spatial resolution (FWHM of the PSF). This en- 
sures that we consider oinly velocities where the line emission 
along the slit is spatially unresolved. 



We first consider a reference frame in the plane of the sky 
which is centered on the center of PA1 slit and which has the 
X axis along the North direction. For a given velocity bin v, the 
position of the light centroid on the sky plane is [x(vf), y(v,-)]; the 
expected photocenter position along a given slit (i.e. a point of 
the spectroastrometric curve) is simply the position [x(v,), y(vf)] 
projected along the slit direction that is: 

PtoM) = Wv,-) - xf]cos{G slil ) + \y(Vi) - yf]sen{8 sW ) (B.l) 

where and y s I lt are the coordinates of the slit centre and 
6 s ii t is the position angle of the slit referred to the X axis (i.e. the 
North direction as in the usual definition of PA). 

The spectroastrometric curve provides the measured centroid 
position along a given slit P'^iyi). Thus, in order to determine 
the free parameters x(v,),.y(vj), (6 s u t is known), one can 

minimize the following^ 2 : 
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(B.2) 



It is worth noticing that the relevant quantities are not the ab- 
solute positions of the slit centers, but the relative ones. Indeed 
we have chosen the center of the reference frame coincident with 
the center of the PA1 slit. The number of free parameters is then 
6 (2 slit centers and the photocenter position on the sky, i.e. 6 
coordinates) compared to 3 data points (the photocenter posi- 
tions from the 3 slits). However, the problem is not undefined 
since many velocity bins are available and the slit center posi- 
tions must be the same for all velocity bins. Therefore, for each 
set of slit center positions x s J: lt ,y s J: lt (slit = 1,2,3) we minimize 
separately all^ 2 . The best xi lt ,y s J: lt values are then those which 
the minimize 

X 1 = 2>?(*o*3tf a ) (B.3) 



if N points are available from each spectroastrometric curve, 
the number of degrees of freedom is then d.o.f. - 3N -4-2N = 



